We are IntechOpen, 
the world’s leading publisher of 


Open Access books 
Built by scientists, for scientists 


5.300 130,000 150M 


ailable International authors and editors Downloads 


Our author among the 


154 TOP 1% 12.2% 


Countries delivered to most cited s Contributors from top 500 universities 


Selection of our books indexed in the Book Citation Index 
in Web of Science™ Core Collection (BKCI) 


Interested in publishing with us? 
Contact book.department@intechopen.com 


Numbers displayed above are based on latest data collected. 
F information visit www.intechopen.com 


ay 


21 


Advanced Numerical Simulation 
of Organic Light-emitting Devices 


Beat Ruhstaller!, Evelyne Knapp”, Benjamin Perucco?, Nils Reinke’, 
Daniele Rezzonico? and Felix Miiller® 

123,47 urich University of Applied Sciences, Institute of Computational 

Physics, 8401 Winterthur 

°6Fluxim AG, 8835 Feusisberg 

Switzerland 


1. Introduction 


Organic light-emitting devices (OLEDs) are novel and efficient light sources that consist of a 
sequence of layers that fulfill distinct electronic and optical tasks. Given the variety of organic 
semiconductor materials available, the improvement of the device performance is a tedious 
and demanding task that often involves numerous experimental optimizations of layer 
materials, thicknesses and sequence. Comprehensive numerical device models have recently 
been developed that master both the physical complexity of the underlying optoelectronic 
processes as well as the numerically challenging system of equations. In this chapter we 
introduce such a comprehensive optical and electronic device model that is able to describe the 
device performance of OLEDs. We first introduce the key device model equations for charge 
and exciton transport as well as the dipole emission model for describing the out-coupling 
of light. Then we present a series of simulation results that are of practical interest when 
studying and optimizing OLEDs. These include the calculation of current-voltage curves, 
current transients signals, time-of-flight current transients and impedance spectroscopy data. 
We show that the physical model can be combined with a nonlinear least-square fitting 
algorithm for extracting transport parameters from measurements. In terms of optical device 
characteristics, quantitative outcoupling mode contributions and angular characteristics are 
presented as well as results of an emission zone extraction method. 

Electrical characterization of devices and materials is essential and helps to elucidate 
the underlying, physical models of charge carrier transport in disordered, organic 
semiconductors. Besides the commonly used current-voltage curves, dark-injection, 
electroluminescence and time-of-flight transient measurements as well as impedance 
spectroscopy offer other ways to validate models for organic LEDs and extract model 
parameters. By means of a one-dimensional numerical OLED model we are able to simulate 
these different measurement techniques. Here, we present numerical methods in the physical 
and numerical framework of reference (Knapp et al, 2010) and solve directly for the steady- 
and transient state. Further, we conduct a numerical small signal-analysis for OLEDs. The 
underlying model solves the drift-diffusion equations in a coupled manner for disordered, 
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organic semiconductors. The disordered nature of organic semiconductors affects the density 
of state, the mobility model, the Einstein diffusion relation as well as charge injection. 
These novel physical model ingredients constitute a second generation OLED model and are 
implemented in the simulator SETFOS (Fluxim AG, 2010). It is expected that the second 
generation OLED model will impact the way OLED characteristics and performance are 
quantitatively described. 


2. Description of the device model 


2.1 Charge drift-diffusion model 

To describe the main features of charge transport in organic LEDs four processes have to be 
considered as illustrated in the schematic energy level diagram in Fig. 1. In a first step, charge 
carriers have to be injected into the organic material (1), secondly they will be transported (2) 
until they recombine to an exciton (3). Then the excitons decay radiatively or non-radiatively 
(4). In the following we will first look at the transport process (2). For the description of 


@ 


ges 


Anode Hole EML Electron Cathode 
transport transport 
layer layer 
(HTL) (ETL) 


Fig. 1. Main processes in OLED operation: 1) Injection, 2) Transport, 3) Formation of exctions, 
4) Radiative decay 


charge transport in OLEDs the general semiconductor drift-diffusion equations for electrons 
and holes are valid. In Poisson’s equation 


V- (eVip) = e(n +n — p = pr), (1) 


the electrical potential yp is related to the mobile electron and hole densities n and p and the 
trapped electron and hole densities n; and p; where e is the elementary charge and e the 
product of the vacuum permittivity €o and the relative permittivity €; of the organic material. 
The current equations for electrons and holes read 


Jn = —enpnVyp+eDnVn, 


2 
Jp = —epupV y} —eDpVp 2) 
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where Hn,p denotes the mobility and Dy,» the diffusion coefficient for electrons and holes. Only 
mobile charges contribute to the current. The conservation of charges leads to the continuity 
equations for electrons and holes 


rae me 6 
F= IVJ, Rap- i 


where R denotes the bimolecular recombination rate given by Langevin and t the time 
(Langevin, 1903). These equations take charge migration and recombination into account. The 
trapped electron (n+) and hole (p+) charge carriers obey the rate equations for an energetically 
sharp trap levels as shown in Fig. 2 


ont — ron(Ni — nt) — fenit, 


Opt 


(4) 
ar = teP( Ni — pt) — rept. 


where re denotes the escape rate, r¢ the capture rate and N; the trap density. Note, that more 
general trap distributions can be introduced that are described by an exponential or Gaussian 
density of trap states (Fluxim AG, 2010; Knapp et al, 2010). 


Fig. 2. Gaussian distributions of density of states and trap levels for trapped charges. 


As opposed to inorganic semiconductors the density of states for organic semiconductors is 
described by a Gaussian shape since transport is assumed to occur via a hopping process 
between uncorrelated sites. Thus, polymers and small molecules have broadened energy 
levels of their highest occupied molecular orbital (HOMO) and lowest unoccupied molecular 
orbital (LUMO) as shown in Fig. 2 and are described in the following way 


2 
No _ (E-Eo 
NGauss (E) = Jina exp | ( Jo ) | (5) 


where No denotes the site density, 7 the width of Gaussian and Eo the reference energy 
level. In the extended Gaussian disorder model the Gaussian density of states affects charge 
diffusion. Tessler pointed out that the use of the generalized instead of the classical Einstein 
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relation is appropriate (Tessler et al., 2002). The generalized Einstein diffusion coefficient is 
now determined by 


D = HT, p,F)gs(p, T), 6) 
where the enhancement function g3 reads 
(PT) = et 7) 
Eş 
Using the expression 
p(Ef) = m DOS(E) f(E, Ep)dE (8) 


for the density and inserting the Fermi-Dirac distribution and the Gaussian DOS we obtain 


f Noauss(E) — rmn E 
D kT” 1+exp ( ET ) 


E=E 
l q a NGauss (E) exp (at 


(9) 


We will now turn to the charge mobility model. A mobility model that has been applied for 
quite some time now is the Poole-Frenkel mobility which is field-dependent and reads 


u = po exp (vv T) (10) 


where ig is the zero-field mobility and y is the field-dependence parameter. However, it has 
been shown by Bassler with the aid of Monte Carlo simulations that the energetic disorder in 
organic semiconductors influences the charge mobility (Pautmeier et al., 1990). Experiments 
have shown that the mobility in hole-only devices can differ up to three orders of magnitude 
between OLED and OFET device configurations with the same organic semiconductor. An 
explanation for this difference is a strong dependence of the mobility on the charge density 
(Tanase et al., 2003). Vissenberg and Matters developed a mobility model that considers 
such a density-dependent effect (Vissenberg & Matters, 1998). Using a 3D master equation 
approach to simulate the hopping transport in disordered semiconductors a dependence on 
the temperature, field and density was determined. Pasveer’s model is therefore dependent 
on the temperature, field as well as the density and accounts for the disorder in the material 
(Pasveer et al, 2005). The extended Gaussian disorder model (EGDM) is an extension of the 
Pasveer model by additionally considering diffusion effects. In the EGDM the mobility may 
be expressed as a product of a density-dependent and field-dependent factor according to van 
Mensfoort as (Mensfoort et al., 2008a) 


u(T, p, F) = wo(T)g1(p,T)g2(F,T), (11) 


with the enhancement functions gı(p,T) and g(F,T) accounting for the influence of the 
charge density, the electric field and the energetic disorder. The functions g1(p, T) and go(F, T) 
are nonlinear and strongly increase with higher disorder parameter (Mensfoort et al., 2008b). 
In the remainder of this text, we will use either the traditional mobility model (Equation 10) 
or the advanced density-dependent model (Equation 11). 
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2.2 Exciton generation, transport and decay 

This section describes the generation, transport and decay of excitons in organic 
semiconductor layers of OLEDs and becomes relevant in bipolar devices. While electron-hole 
recombination may be locally confined, the resulting excitons can migrate with a characteristic 
diffusion constant before decaying radiatively or transferring their energy to other exciton 
species in the same layer. Thus, the rate equation for exciton species S; contains a generation, 
a diffusion, transfer and (non-)radiative terms: 


dS;(z) 
dt 


= G;R(z) + VJs, (2) = (krad; (Z) + Knonrad; ) i S;(z) Ti Kannihilation; à Si(z)? F (12) 


exc 


j= 


In Eq. 12 G; is the generation efficiency prefactor, k,,q,(z) the position-dependent radiative 
decay rate, knonrad; the non-radiative decay rate and Kkannihilation; an annihilation rate. For 
a comprehensive fully-coupled opto-electronic OLED model, the radiative decay rate is 
calculated with the dipole emission model (see section 2.3 below) and ref. (Fluxim AG, 2010). 
In a classical view the generation efficiency for singlet excitons is G; = 0.25 while it would 
be G; = 0.75 for triplet excitons. The exciton energy transfer term —k;j - S; shall be used to 
model energy transfer from exciton S; to exciton species Sj. Such exciton energy transfer rates 
can be taken from the literature. Common energy transfer mechanisms are Förster and Dexter 
transfer (Dexter, 1953; Förster, 1948). Forster transfer is a long-range dipole-dipole interaction 
while Dexter transfer is a short-range electron-exchange interaction. Both mechanisms rely 
on the spectral overlap of the emission spectrum of the donating exciton species with the 
absorption spectrum of the accepting exciton species. Exciton transport happens by diffusion 
and is therefore mathematically described by 


Js; = D:V Si (13) 


where D, denotes the exciton diffusion constant. 


2.3 Optical dipole emission model 

In this section we present a model for the dipole emission describing the radiative decay 
of excitons in planar multi-layer systems. This model accounts for the influence of layer 
interfaces on the radiative decay rate of fluorescent molecules and is based on a classical 
theory for the emission of an oscillating electrical dipole. The basic concept of this theory 
was derived by Sommerfeld to analyze electromagnetic radiaton close to a conductive surface 
and extended by Chance et al. (Chance et al, 1978) to model fluorescence emission in optical 
cavities. 

The emission characteristic of a single molecule is similar to that of a Hertzian dipole. 
However, in an infinite medium consisting of isotropically oriented molecules, the emission 
from the multitude of molecules is homogeneous in all directions. If, on the contrary, the 
molecules are embedded in a thin layer, interference effects of electromagnetic waves reflected 
at layer interfaces dominate the emission characteristics. 

The dipole emission model allows to quantitatively analyze the individual energy loss 
channels of excited molecules. These energy loss channels include radiative coupling to the 
far-field, modes guided inside the multilayer structure as well as absorption losses by the 
surrounding media. The model is based on the following assumptions: 
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(1) Each medium of the multilayer system is assumed to be isotropic and only characterized 
by its dielectric permittivity e€. 


(2) The medium embedding the emissive dipole is considered to be loss-free (i.e. real 
dielectric permittivity e€) 


(3) The interface planes are infinitely extended. 


(4) The wavelength A of the radiation as well as the distance between the emitter and the 
individual interfaces is large with respect to the spatial extensions of the dipole (point 
dipole approximation). 


(5) The dipole is at least at A /50 distance away from the regions with a complex permittivity 
€. 


(6) The materials of the multilayer stack are considered to be non-magnetic. 


Despite these restrictions, this model is valid for a wide range of physical situations. 
The model considers the molecules as driven damped harmonic dipole oscillators. These 
oscillators are damped by the emission process and driven by electric field components 
reflected at dielectric interfaces. An illustration of an emitting dipole near an interface is given 
in Fig. 3. The dynamics of the oscillating dipole moment # is described by the differential 


Fig. 3. Illustration of an emissive dipole near a reflecting interface at distance d. 


equation 

d2 dp 2 e2 > 

“RR “ike f— —F | 14 

qa +t bog, to Pp = ER w) (14) 
with the oscillator frequency w, the elementary charge e, the effective mass m of the dipole, 
and Er(w) the electric field at the dipole position. The radiative power of a dipole located 
inside an isotropic, infinite medium with a refractive index of ne is 


2wk3 
qobo = pL hy 


127€9 ~’ 


where go is the intrinsic quantum efficiency of the dipole and kg = 27t/Apo. The intrinsic power 
of a dipole consists of a radiative and a non-radiative component. 

The dynamics of a dipole located in an optical multilayer system is affected by the presence 
of the reflected field Er and from its radiative character as 


b = gb + (1 — q)b, (15) 
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where q is the modified quantum efficiency. The second term in Equation 15 refers to the 
non-radiative processes within the optical structure. These are not affected by the optical 
environment. Therefore we summarize the impact of the optical environment on the radiative 
part and write 

b = gobo - F + (1—qo)bo (16) 
where F describes the influence of the electromagnetic field components on the oscillating 
dipole. F can be expressed by means of a generalized parametrization of the emission angle 
u = sin ĝe as 


Ee f fooau. (17) 


The total lifetime of a dipole in the optical system is consequently related to its intrinsic 
lifetime To by 
TO b 


ooh ae (18) 


Please note, that the above equations are modified for distinct dipole orientations and 
polarizations of the emitted light. 


2.3.1 Mode Analysis 

As seen above, the total dissipated power from the dipole is proportional to F and therefore, a 
stepwise integration of Equation 17 based on the index limits of total internal reflection allows 
to determine the fractional radiated power that can be assigned to out-coupled, guided (in the 
thin films), and evanescent modes. The integration limits are given by: 


u = 0 snil Ty radiative modes (19) 
= Wis Bens | guided modes 


=1...0 evanescent modes, 


nt and ne being the refractive indices of the top and the emissive medium, respectively. The 

multilayer surrounding of the emitting dipole does not only affect its lifetime T = Tọ - bg/b, 

but also its radiative quantum efficiency. From Equation 15 and 16, the modified quantum 
efficiency q is therefore 

qo: F 

= =. 20 

Dewey (20) 


3. Simulation results 


In this section, we present a number of simulation results that i) highlight the impact 
of energetic disorder on steady-state and transient currents as well as impedance spectra 
(subsections 3.1 to 3.3), ii) exemplify the fully-coupled electronic-optical simulation of 
multilayer OLEDs (subsection 3.4) and iii) demonstrate the use of time-of-flight transients 
under the influence of traps (subsection 3.5). 


3.1 Steady-State Current-Voltage Curves 

In the following we present the simulation results of a symmetric hole-only device under 
three different operating conditions. Thereby we consider the disordered nature of organic 
semiconductors that affects the density of states, the mobility model, the Einstein relation 
as well as charge injection. In this device the highest occupied molecular orbital (HOMO) 
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is aligned with the workfunctions of the anode and cathode. Therefore, we obtain a 
well-injecting anode and cathode. The high barrier for injection of electrons from the cathode 
to the lowest unoccupied molecular orbital (LUMO) level essentially leads to a hole-only 
device. 

The high charge densities at the anode and cathode (p(0) and p(L)) are set to half of the 
density of chargeable sites No, i.e. the Gaussian DOS is half-filled at either side due to the 
energy level alignment at the electrodes. 


Parameter Value Units 


No 8.5- 1076 m~’ 

o 0.15 eV 

Ho 10710 m2Vy~!571 
Voi 0 V 

p(0) 0.5No m~ 

p(L) 0.5No m? 


Table 1. Simulation parameter set of the hole-only device for the three operating conditions. 


First, we look at the best known electrical characterization method: the current-voltage curve. 
Therefore the steady-state for different voltages is calculated and we obtain current-voltage 
characteristics as shown in Fig. 4. We compare the steady-state results for the EGDM and the 
constant mobility case with the analytical solution which is given by 


9 y? 
J(V) = g EEO HO T3 (21) 


where L denotes the length of the device, V the applied voltage and pig the mobility (Mott 
& Gurney, 1938). All of the analytical solutions in the following are obtained by neglecting 
the diffusion and only taking the drift in Equation 2 into account. We observe the effect of 
diffusion especially at low voltages where the current is increased. In the regime of high 
bias the analytical solution and the constant mobility case merge whereas the EGDM solution 
moves away due to the field- and density-dependence of the EGDM and the slope of the 
EGDM curve is increased. Generally, the current density of the EGDM is higher due to 
the enhancement functions 91(p,T), &2(F,T) and g3(p,T). The EGDM has recently been 
successfully applied to current-voltage curves by van Mensfoort et al. (Mensfoort et al., 
2008a). A comprehensive analysis of the effects of the EGDM and its components for the 
steady-state is conducted in reference (Mensfoort et al., 2008b). 


3.2 Dark-injection transients 

We now turn to the analysis of the time-dependent response to a step voltage. We perform 
transient simulations for the same device as above shown in Fig. 5. The dark-injection 
transients are simulated from 2 to 12 V in steps of 2 V. For the initial rise of the transient 
response the analytical solution is described by as 


_ 2Leeg 1 
y (2tscr — t)?’ 


J(t) (22) 


where the space-charge-free transit time is defined as tscr = 4p (Helfrich & Mark, 1962). In 
the SCLC at the transit time Tt = 0.786Tscr the maximum current is obtained and reaches 
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= = = const 
analytical 


current density J [A/m] 


107 107 10 10 
voltage [V] 
Fig. 4. Current-voltage curve for a hole-only device for a constant mobility and diffusion 
coefficient, the EGDM and the analytic solution (Mott-Gurney). 


a value of J(t%) = 1.21Je where Je stands for the stationary current from Mott-Gurney’s law 
(Equation 21). The value of J(t) = 1.21J, is marked in Fig. 5 with a star. Note that the slope 
of the line connecting these stars is minus two which can be explained with Equation 21. 
The relatively low charge mobility of organic semiconductors leads to space-charge limited 
transport. In this regime the transient response exhibits a cusp that relates to the transit time 
and thus the charge mobility. This has been the motivation for carrying out dark-injection 
transient experiments. 

We notice that the peak position for the constant mobility case coincides well with the 
analytical solution, especially at high bias where the drift current is the dominant part of the 
total current. At low voltage the agreement becomes worse due to the diffusion current. 


e 1.21 
= = =const 
analytical 


current density J [A/m?] 
3 


time [s] 


Fig. 5. Dark-injection transients for a hole-only device for a constant mobility and diffusion 
coefficient, the analytical solution and the generalized Einstein solution. The transients were 
simulated from 2 to 12 V. 
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Fig. 6. Dark-injection transients for a hole-only device for the EGDM at different voltages. 


If we also consider the enhancement functions for the EGDM as shown in Fig. 6, the peak 
positions shift to shorter times due to the enhanced mobility and the current density is vastly 
increased. In the EGDM case no formulas are available that would relate the peak position in 
the current with EGDM mobility parameters. Therefore, we have to rely on numerical fitting 
of the model parameters. 


3.3 Impedance Analysis 
As a third step, we investigate the response of the device to a voltage modulation. We obtain 
therefrom the impedance Z(w) which is related to the admittance Y(w) as 


Y(w) = ——~. (23) 


The admittance is decomposed into an equivalent parallel conductance and capacitance as 
follows 
Y(w) = G(w) + iwC(w). (24) 


We analyze the normalized capacitance at 8V in Fig. 7. We varied the energetic disorder 
parameter 7 which is the main EGDM model parameter. We compare our numerical 
simulation results with the analytical solution for the drift-only case. For the small signal 
analysis the analytical solution for the impedance Z(w) is given by 


6 


A= goliwTo)’ 


[1 — iwTo + Lliw) — exp (—iw7)| (25) 


where To = T and gy = 2e (Shao & Wright, 1961). We also calculated the solution 
for the constant mobility case. The drift-only model provides a quite accurate solution for 
the constant mobility case at 8 V. At lower voltages however, the drift-only model differs 
significantly from the constant mobility model. We observe that the oscillations at high 
frequencies are damped out when diffusion is present. Also due to diffusion, the capacitance 
at low frequencies is increased from the typical value IC seoj of the drift-only solution, even 
more in case of the EGDM. In the EGDM the peak position is additionally shifted to higher 
frequencies. 
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Fig. 7. Frequency dependence of the normalized capacitance at 6 V in a hole-only device for a 
constant mobility and diffusion coefficient, the EGDM and the analytical solution. 


3.4 Fully-coupled electrical-optical OLED simulation 

State-of-the art OLEDs are composed of several functional layers whose thickness obeys 
some optical or electronic criteria and which may be further complicated by dye doping, 
charge doping etc. The efficient harvesting of excitation energy by use of fluorescent and 
phosphorescent dye dopants remains a challenge. Here we attempt to give a simulation 
example of dual emitters in two layers EML1 and EML2 of a hypothetical 4-layer device 
structure. The rate equation for excitons (see Equation 12) is used here to model the energy 
transfer among different exciton species in a multilayer small-molecule OLED as depicted in 
Fig. 8. 
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Fig. 8. Energy level diagram of a 4-layer OLED with two emissive layers EML1 and EML2 
with two exciton species each. 


The middle two layers were assumed to have emissive excitons in the host and guest material. 
The mobilities were assumed as illustrated in Fig. 10. The resulting charge distribution is 
shown in Fig. 9. For the chosen parameter values, the charges accumulate at either side of the 
EML2/ETL interface, thus also leading to a high electric field at that interface, see Fig. 9. 

While only the excitons of the host materials EML1 and EML2 are excited by recombination 
events, the guest dye dopants receive their energy from the host by transfer. For this energy 
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Fig. 9. Simulated profiles of charge carrier densities (left) and electric fields (right) for the 
hypothetical 4-layer OLED example of Fig. 8. 


transfer we employ the transfer rate kj; in Equation 12) and find that the guest excitons 
equilibriate in a different spatial distribution, due to a different parametrization, see Fig. 10 
(right). The excitons are only modeled in the middle two layers. 
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Fig. 10. Simulated profiles of mobilities (left) and Exciton.1 and Exciton.2 in the middle of the 
OLED structure (right). 


If one is only interested in the steady-state, then the governing transport equations 3 are 
solved with an iterative steady-state solver (Fluxim AG, 2010; Knapp et al, 2010). However, 
if one is interested in the dynamics of charge transport, transient solution algorithms must 
be employed. Experimentally, the transient electroluminescence and transient current signals 
have revealed the dynamics of charge migration into the device (Ruhstaller et al., 2003; 2001). 
In Fig. 11 one can observe a monotonic increase of the recombination current, which testifies 
that the charges have first to penetrate the device before they can recombine. For the given 
hypothetical device, the charges are all consumed in the device, i.e. the recombination current 
is as large as the device current and therefore the current balance is equal to unity. 
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Fig. 11. Simulated device current and recombination current during turn-on. 


This example simulation can be completed with a light-outcoupling calculation which takes 
the exciton distribution as dipole distribution. In this step one has to select which excitons 
are assumed to decay radiatively and thus contribute to light emission. For instance, we 
consider the emissive guest excitons in the emission layer EML1 and perform a calculation 
of the dissipated power versus the wavelength and the normalized in-plane wave vector u in 
Fig. 12 (left). The various emission modes can easily be identified qualitatively by looking at 
this plot. The guided modes are visible in the u range between 0.6 and 1. This power seems to 
exceed the power in the propagating outcoupled mode (u < 0.6) while the evanescent modes 
are displayed at values of u larger than 1. In addition in the left plot of Fig. 12 we report an 
analysis of emissive modes that takes care of the integration of the dissipated power of each 
mode in the relevant u ranges. The critical in-plane wavevector for light outcoupling into air 
is around 0.6 and the corresponding relative mode contribution is not more than 10%, in this 
example. Outcoupling efficiency numbers around 20% are typical values for OLEDs without 
any outcoupling enhancement techniques. 


3.5 Numerical analysis of time-of-flight measurements 

An exhaustive study of transport mechanisms in organic compounds shall comprehend the 
characterization of the charge dynamics. Time-of-flight (TOF) experiments are among the 
most common techniques to measure the mobility of the charges. Nevertheless, only a 
thorough analysis of the measured transient curves can give insight into the behavior of the 
charges, including the impact of trapped states. Current transients in TOF measurements are 
generally classified into dispersive and non-dispersive ones, the former being characterized 
by a monotonic decay of the measured current, while the latter show an initial plateau 
followed by a clear and sudden drop in the current. 

In the following, we will investigate TOF transients by numerical simulations for distinct 
kinds of samples in order to allocate the specific characteristics of electron trap states. Here, we 
reproduce the TOF experiment presented by Malliaras et al. for Alq3 devices that underwent 
different compound purification steps (Malliaras et al., 2001). The samples of our virtual 
experiment consisted of 8-~m Alq3 films sandwiched between a transparent ITO layer and 
an aluminium (Al) electrode. The samples were reverse biased at a voltage of —1000 V, the 
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Fig. 12. Dipole emission calculation for the guest excitons in EML1. The dissipated power 
(left) and the mode analysis (right) are shown. 


ITO being the negative electrode. A 10 ps laser pulse at 337 nm wavelength with energy of 1pJ 
per square meter was used to generate a sheet of free carriers at the ITO-Alq3-interface. The 
light intensity was kept sufficiently low to avoid space charge effects. Assuming HOMO and 
LUMO levels for the organic compound to be 5.6 and 3.0 eV, respectively, the two electrodes 
were acting as non-injecting layers because of their workfunctions of 5.0 eV for ITO and 4.1 eV 
for Al. For completeness, the (vanishing) injection was considered as thermionic and modeled 
after Scott and Malliaras (Scott & Malliaras, 1999). 

Despite the energetic disorder of the transport states in organic compounds, the analysis was 
performed assuming Poole-Frenkel, field-dependent mobility, see Equation 10. The study of 
the impact of the energetic disorder on the current transient is beyond the goal of this text 
and it is known that the disorder does not affect the dispersive characteristics of TOF-related 
current signals (Bassler, 1993). Trapped charge carriers on sharp energy levels were described 
by Staudigel et al. (Staudigel et al, 1999) by a rate equation for trapped electrons n; that reads 


d kT Etn 
OD) a aene (26) 
dt edm 


In Equation (26), No is the total number of available sites per unit volume, N; is the trap 
concentration, E is the (local) electric field while E;, is the trap energy depth. Further, dm = 
N; 1/3 is the mean trap distance and o; = 5d7, is the effective capture cross section while n 
is the concentration of free carriers. Note, that Equation 26 is connected to Equation 4 by the 


relation re = Eo; and re = cNoe` T., The electrical parameters for the modeled compound 
are summarized in Table 2. All calculations were performed with SETFOS (Fluxim AG, 2010). 
In Figure 13 we present the results for calculations considering no traps, shallow and deep 
traps. 

When no traps are considered, a non-dispersive charge behavior is clearly recognized in the 
current transient. Contrarily, when deep traps (typically E;, > 0.3 eV at room temperature) 
are taken into account, the transient characteristics is evidently dispersive since no definite 
current drop is recognizable when the free electrons reach the aluminium electrode. A more 
careful analysis is needed in case of shallow traps (E;, < 0.25 eV at room temperature). Due 
to the relatively low energy barrier, trapped electrons are quickly released into the conduction 
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n p 
uo 2.9e-9 cm^2/Vs |1e-8 cm^2/Vs 
y 7.3e-3 cm^ 0.5/V^ 0.5 0 
LUMO/HOMO 3.0 eV 5.6 eV 
A — /0.212/0.33 eV = 


Table 2. Assumed electrical parameters for Alq3 to reproduce the measurements presented in 
reference (Malliaras et al., 2001). Three different electron trapping depths, E;, were 
considered. 


band. The resulting current transient characteristics shows an evident drop as the charges 
reach the aluminium electrode. Such “quasi-non-dispersive" behavior is often confused with 
trap-free characteristics (Malliaras et al., 2001), whereas trap states are still present but do not 
significantly damp the dynamics of the charges. 


No traps 
—-—-- Shallow traps, E;,,=0.212 eV 
— - — Deep traps, E;,=0.33 eV 


Current [mA/cm?] 


80 100 

Time [ps] 
Fig. 13. Current transient characteristics of TOF measurements for samples with no traps (red 
solid line), shallow traps (blue dashed line) and deep traps (green dot-dashed line). As 
expected, trap-free compounds show a non-dispersive charge transport and charges in 
materials with deep traps exhibit a dispersive behavior. Weakly trapped charges on the 
contrary would indicate a non-dispersive behavior which is often mistakenly confused with 
the absence of traps (blue dashed line). 


The behavior of the electrons in presence of deep or shallow traps can be analyzed looking at 
the charge density profiles for different times. Immediately after the light pulse that generates 
the charges a package of free electrons moves away from the ITO layer toward the positively 
charged aluminium electrode. 

In Figure 14, left, we show the profiles of the charges in presence of deep traps, 10 pis after 
the light pulse. We recognize a symmetric, free electron packet partially overlapping the 
profile of the trapped electrons. When the free charges reach the aluminium electrode after 
60 ps, their number is strongly reduced due to the many trapped states occupied during 
the flight from one electrode to the other (Figure 14, right). For the trapped electrons, it is 
energetically unfavorable to release back into a conductive state. They need many hundreds 
of microseconds to be de-trapped. 
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Fig. 14. Deep traps. Left: At t = 10 us, the sharp packet of free electrons (red solid line) leaves 
the ITO electrode at position 0 um and moves toward the aluminium electrode placed at 
position 8000 um. A certain amount of charges remains trapped close to the ITO electrode 
(blue dashed line). Right: After 60 us, those few free carriers that have not been trapped 
reach the aluminium electrode. Notice that trapped charges were “left behind" by the 
traveling free charge packet and did not move from their original position. Few trapped 
charges that begin to release back into the conduction band are recognizable through the 
slight increase of the free charges concentration around 2000 pm. 
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Fig. 15. Shallow traps. Left: Already at t = 10 pis, the free charge packet (red solid line) 
appears asymmetric due to the release of the weakly trapped charges (blue dashed line) 
around 1000 um. Right: After 60 us, a significant amount of free charges reaches the 
aluminium electrode. The fast kinetics of the trapped states makes them release back into the 
conduction band fast enough to follow the movement of the free electrons. 
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In Figure 15, left, we show the profiles of the charges in presence of shallow traps, 10 ps after 
the light pulse. The profile of the free electrons appears asymmetric because of the on-going 
release process of the temporarily trapped electrons. The kinetics of the trapped charges is so 
fast that the trapped packet follows the free electron packet in its travel toward the positive 
electrode. In Figure 15, right, we show the charge profiles at the moment they reach the 
aluminium electrode after 60 us: the broad free electron packet overlaps the trapped charge 
packet so that both types appear to be extracted about at the same time. This explains the 
non-dispersive-like characteristics of the current transient in case of trapped states that are 
less than 25 meV below the LUMO level (blue dashed curve in Figure 13). 

Despite the lack of precise information on the experimental procedure, the calculated current 
transients reproduce both qualitatively and quantitatively the measured curves presented by 
Malliaras et al.(Malliaras et al., 2001). Through our comprehensive numerical investigation, 
we could show that a qualitative analysis of the shape of the current transient is not enough to 
characterize trap-free samples. As a matter of facts, weakly trapped charges are able to follow 
the “flight" of free electrons through the device. Underestimating the presence of traps might 
lead a.o. to wrong statements on the lifetime of the device. 


4. Parameter extraction 


The following section introduces some aspects of parameter extraction in organic 
light-emitting devices (OLEDs). Two main applications are discussed. The first application 
is the extraction of the emission profile in a light-emitting layer of an OLED by a numerical 
method. The emission profile describes mathematically how the excitons are distributed 
inside the emissive layer. The numerical method presented is a linear least-square fit 
algorithm. The theory of numerical algorithms to extract the emission profile and some 
applications are also presented in more detail in Perucco et al. (Perucco et al, 2010). The second 
application is the extraction of EGDM parameters from multiple measured current-voltage 
curves by a nonlinear least-square algorithm. 


4.1 Extraction of emission profiles in OLEDs 

The objective of this section is to present and test a numerical fitting algorithm for the 
extraction of the emission profile and intrinsic source spectrum. The fitting algorithm is 
evaluated by adequate examples and validated on the basis of consistency checks. This 
is achieved by an optical model, where a transfer-matrix theory approach for multi-layer 
systems is used in combination with a dipole emission model. The optical model is 
implemented in the semiconducting emissive thin film optics simulator (GETFOS) (Fluxim 
AG, 2010). With SETFOS, the emission spectrum of an OLED based on an assumed emission 
profile and a known source spectrum is generated. The fitting method is then applied to the 
calculated emission spectra in order to estimate the emission profile and source spectrum. 
The comparison between the obtained and assumed emission profile and source spectrum 
is an indication of how successfully the inverse problem can be solved. Sections 4.1.1 deals 
with the mathematical derivation of this numerical fitting algorithm. In Section 4.1.2, the 
applications or consistency checks are presented. 


4.1.1 Theory 

The theoretical background of the fitting method is introduced in this section. The method 
is linear in terms how the measured emission spectrum is related to the unknown emission 
profile. For simplicity the mathematical formulation for the least-square problem is derived 
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for only one emitter. The emitter is characterized by the emission profile. For the moment, it 
is also assumed that the emission spectrum is measured in normal direction and therefore 
the light is unpolarized. In Section 4.1.1.1, this approach is extended further to multiple 
emitters described by several emission profiles and emission spectra measured for several 
angles 0. The extracted emission profile P.(d;) is discretized at N relative positions 6; = dj/L 
in the light-emitting layer, where dj is expressed as an absolute position and L is the width 
of the layer. The emission spectrum is divided into M wavelengths A; (i = 1...M). The fitted 
emission spectrum I(A;) can be written as 


N 
Ip(Ai) = Yo Ie (Aj, 6)) - Pe(5j), (27) 


where I-(Aj;,6;) is the emission intensity for the wavelength A; and assuming a discrete 
emission profile (dirac function) at the relative position 4; in the layer. The emission intensity 
is given by 

To(Aj, 8j) = [(Aj, 6;) - Si), (28) 
where I(A;, ôj) is the emission intensity for emissive dipoles with spectrally constant intensity. 
S(A;) is the source spectrum. Between the measured emission spectrum Im(A;) and fitted 
emission spectrum I 6 (Ai), a residuum can be defined and written as 


(Aj) = If (Aj) — Im(Ai) . (29) 


Equation 29 can be interpreted as a linear least-square problem, written as a system of linear 
equations 


N 
r1(Ai) = Yo Ie(Ai 8j) < Peð) — Im(Ai) - (30) 
jal 


The system of equations is normally overdetermined (i.e. M > N) and thus is ill-posed. In 
matrix notation, the problem can be formulated as r4 = A- xı — bı. The matrix A has the 
following structure 


I.(Am, 61) I-(Am, dy) eins I-(Am, On) 


bı is a vector containing the measured emission spectrum Im(A;) and the vector xı 
corresponds to the a priori unknown emission profile P.(d;). The term linear refers to the 
linear combination between the matrix A and the vector x; of unknown weights. In every 
column of the matrix A, an emission spectrum is calculated for a dirac shaped emission profile 
at the position 6;. The emission profile P,(0;) at the relative position 4; is the weight of the 
corresponding spectrum, respectively the column. The mathematical task is to minimize the 
length of the vector || rı ||. 


4.1.1.1 Extracting multiple emission profiles 


The most general case of the emission spectrum is determined by the emission profile of 
multiple emitters P.(5") and emission angles 6;. Given is the emission spectrum measured 


at O different angles (J = 1...0) and the OLED consists of Q different emitters (k = 1...Q) in 
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the same or in separate layers. The relation stated in Equation 27, combined with the definition 
of the residuum in Equation 30, can be extended to 


P (Ai, 01) = yer (Ai, OF, 01) + P(E) — In? (Aj, 81) - (32) 


EP (A; of ,0,) stands for the s-polarized or p-polarized emission intensity at the wavelength 
A;. We assume a dirac shaped emission profile at the relative position of for emitter k and an 
emission angle of 6). Pe (ð; ) is the emission profile at relative position 6; for emitter k. Equation 


32 represents a system of linear equations r3” = A®? - x. — by”, where the matrix A®? contains 
the s-polarized and p-polarized emission spectra, the vector x2 contains the information of 
several emission profiles and the vector bP represents the measured emission spectrum. The 


mathematical task is again to minimize the length of the vector || r3” ||. 


4.1.1.2 Extracting the intrinsic source spectrum 


In the case of a single emitter, van Mensfoort et al. (Mensfoort et al., 2010) presented a 
method to extract the source spectrum of the light-emitting material. The source spectrum 
can be obtained by replacing the emission intensity I’? (Aj, ot, 01) by the emission intensity for 


emissive dipoles with spectrally constant intensity [°? (Aj, ot, 01) in Equation 32. This method 
is employed and evaluated in Section 4.1.2.2. 


4.1.2 Applications 

In this section, the reliability and limitation of the linear fitting method is addressed after it 
was mathematically deduced and described in Section 4.1.1. A given intrinsic source spectrum 
from a light-emitting material is assumed, together with an emission profile, stating where the 
dipoles are located in the device. The effects of quenching are disregarded in the presented 
applications below. First, quenching would likely limit the amount of dipoles close to the 
electrodes as the lifetime is very short. And secondly, light emitted from the dipoles is also 
captured in evanescent modes and therefore, does not couple out into air. Finally, the emission 
spectrum is generated by an optical dipole model described by Novotny (Novotny, 1997) and 
implemented in the simulator SETFOS (Fluxim AG, 2010). The calculated emission spectrum 
is used to solve the least-square problem in Equation 32. This allows the extraction of both, 
source spectrum and emission profile. The comparison of the extracted and assumed emission 
profile reveals the reliability of the presented algorithm. Throughout this text, an open cavity 
is used for the consistency checks. But the method here may also be applied to cavity 
and small-molecule based OLEDs. The OLED investigated here has a broad light-emitting 
polymer (LEP) of 100 nm. Further, the light-emitting layer is embedded between a 80 nm thick 
PEDOT:PSS anode and an aluminum cathode of 100 nm. The device is depicted in Figure 16. 
With respect to an experimental setup, the diameter of the semi-sphere glass lens is at least 
an order of magnitude larger than the diameter of the OLED. In order to achieve an absolute 
quantity of the emission intensity and emission profile, the assumed current density in all 
considered consistency checks is 10 mA/cm?. 


4.1.2.1 Extraction of the emission profile from angularly resolved emission spectra 


As an introductory example, this section shows the application to angularly resolved emission 
intensity spectra. It compares the extracted emission profile from these spectra to an emission 
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Fig. 16. OLED used to perform the parameter extraction tests with the semi-sphere glass lens. 
0 stands for the observation angle. 


profile extracted from an emission intensity spectrum measured at normal angle. The 
assumed emission profile is Gaussian shaped, where the peak is set to 0.3 expressed in terms 
of a relative position in the emission layer. The width of the Gaussian shape is 20 nm. 
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Fig. 17. Comparison between the assumed and extracted emission profiles. The emission 


profiles were extracted using angularly resolved emission spectra and an emission spectra 
measured at normal angle. 


The comparison between the extracted and assumed emission profiles in Figure 17 shows an 
improvement of the extracted emission profile when angularly resolved emission intensity 
spectra are used. The fitted emission intensity spectra match visually perfectly the emission 
spectra serving as a measurement, as seen from Figure 18. 


4.1.2.2 Source spectrum extraction 


This section demonstrates the ability of the least-square algorithm to extract the intrinsic 
source spectrum of a light-emitting material. The same assumptions regarding the parameters 
of the emission profile are made as in Section 4.1.2.1. The extracted emission profile and source 
spectrum can be found on the left, respectively on the right in Figure 19. 

It can be seen from Figure 19 that the source spectrum can be extracted very accurately. 
The emission profile is also well extracted and even the peak position is reproduced well. 
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Fig. 18. Left: Angularly resolved emission intensity spectra serving as a measurement. Right: 
Fitted emission intensity spectra by the linear least-square algorithm. 
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Fig. 19. Left: Comparison between the assumed and extracted emission profile. Right: 
Relation between the assumed and extracted intrinsic source spectrum by the method 
discussed in Section 4.1.1.2. 


Illustrated in Figure 20 is the comparison between the assumed and fitted emission spectra, 
which are visually also in perfect agreement. 


50 "5 50 "5 
80 457, 80 457 
70 40E 70 40 E 
> 60 35% > 60 35% 
D, 50 30 D, 50 30 £ 
® 40 25= © 40 2a. 
D 202 © 20 = 
© 30 @ S30 O 
& 152 < 15 2 
20 2 20 o 
10 = 10 £ 
10 5 = 10 5 £ 
0 0 a 0 0 = 
520 540 560 580 600 620 640 660 520 540 560 580 600 620 640 660 
Wavelength [nm] Wavelength [nm] 


Fig. 20. Left: Angularly resolved emission intensity spectra serving as a measurement to 
extract the intrinsic source spectrum from. Right: Fitted emission intensity spectra by the 
linear least-square algorithm. 
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4.1.2.3 Extracting multiple emission profiles 


Equation 32 explains how multiple emission profiles can be extracted from a measured 
emission intensity spectrum. This section illustrates the application of the method to a 
multi-emitter OLED. In this example, two emission profiles are extracted. The first assumed 
emission profile is Gaussian shaped with a peak at 0.7 and a width of 40 nm. The second 
assumed emission profile is also gaussian shaped, where the peak is at 0.3 and the width is 
20 nm. Figure 21 shows the extracted and assumed emission profiles, as well as the reference 
and fitted emission intensity spectra. 
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Fig. 21. Left: Showing the differences between the emission profiles from a multi-emitter 
OLED. The dotted curves represent the assumed emission profiles, whereas the lines stand 
for the extracted emission profiles. Right: Comparison between the measured and fitted 
emission intensity spectra. 


Despite the fact that the assumed and fitted emission spectra match very well, some 
differences in the emission profiles are visible. Nonetheless, the general trend is explained 
by the extraction. For the second emission profile, both peak and width can be reproduced 
more or less. For the first emission profile, the flat emission profile can be explained as well. 


4.2 Extraction of transport parameters from current-voltage curves 

The following section deals with the application of a nonlinear least-square fitting algorithm 
to extract EGDM parameters from measured current-voltage curves. The nonlinear fitting 
algorithm, as well as the EDGM model is implemented in SETFOS. SETFOS is used to generate 
three hypothetical measured current-voltage curves at temperatures 320 K, 300 K and 280 K. 
All three current-voltage curves are simultaneously fitted for extracting the parameters. The 
device considered is a single-layer, hole-only device where the electrical layer has a thickness 
of 121.5 nm and the build-in voltage is 1.9 V. The energy diagram of the simulation device is 
depicted schematically in Figure 22. The following parameters are of interest: the mobility up, 
the width of the Gaussian DOS gp, the density of chargeable sites No and the workfunction at 
the cathode ®,. Meanwhile, the workfunction at the anode is held constant. The parameters 
represent real EGDM parameters as discussed in van Mensfoort et al. (Mensfoort et al., 2008b). 
The following parameters are assumed: up = 1-107” m*/Vs, op = 0.13 eV, No = 6: 107° 1/m° 
and ®, = 3.2 eV. The mobility up is related to Equation 11 in the following way: vo(T) = 


Lp exp (—0.39(0/ (kT) )?). The left hand side of Figure 23 shows that the nonlinear least-square 
algorithm is capable of extracting all four EGDM parameters as the current-voltage curves at 
the same temperature match each other visually perfectly. 
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Fig. 22. Energy diagram for the simulated device. The workfunction at the anode side is held 
constant at 5.1 eV, while the workfunction at the cathode side ®; is being optimized. The 
HOMO and LUMO levels are 5.2 eV, respectively 2.1 eV. 
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Fig. 23. Left: Comparison between the measured current-voltage curves (points) and the 
fitted current-voltage curves (lines). The extraction of parameters is done using all three 
measured curves at different temperatures at the same time in order to get more accurate 
results. Right: Showing the relative change in the parameters as a function of iterations. The 
black line represents the relative change which must be achieved as a termination criterion. 


The right hand side of Figure 23 illustrates the relative error of the EGDM parameters as a 
function of the optimization iteration. The relative change is defined as (pj41 — pi)/Pi+t, 
where p; is the value of parameter p at iteration i. The black line shows the relative change, 
which must be achieved in order to terminate the optimization process. All four parameters 
obtained after 9 iterations have a smaller relative change than 10710. Depending on the 
parameter, it can be seen from Figure 23 that the relative change in the parameter is initially 
between 1071 and 1. The assumed, extracted and initial parameters for the fit can be seen in 
Table 3. 

Table 3 shows that the nonlinear least-square algorithm is able to accurately extract the 
assumed EGDM parameters. A potential limitation of the nonlinear least-square method is 
that the error landscape should be smooth and the initial parameters should not be too far off 
from the assumed values, otherwise the method might converge into another local minimum. 
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Parameter Assumed] Extracted | Initial 
Up [10-7 m7/Vs]| 1.400 | 1.398 | 3.000 
Oy [eV] 0.130 0.130 | 0.100 
No [1076 1/m3] 6.000 5.984 | 1.000 
®, [eV] 3.200 3.200 | 3.200 


Table 3. Assumed and extracted EGDM mobility model parameters for the test, as well as 
initial parameters for the nonlinear least-square algorithm. 


For such situations and to cross-check the fit result, potentially more time-consuming global 
optimization routines may be used (Fluxim AG, 2010). 

In a post-processing step, the correlation coefficients for the four model parameters can be 
calculated at the best fit point. These coefficients corresponding to the fit result listed in Table 
3 are shown in Table 4. 


Hp Op No P. 
Hp |1.00/0.94/0.39|-0.95 


Tp 1.00) 0.07} -0.84 
No 1.00)-0.55 
dD, 1.00 


Table 4. Correlation coefficients for the EGDM parameters at the best fit point in case of a 
multiple current-voltage curve fit as seen in Figure 23: mobility yp, width of Gaussian DOS 
Jp, density of chargeable sites No and workfunction cathode ®, 


A correlation of 0 indicates no correlation between the two parameters, meaning the 
parameters are independent from each other. Whereas a correlation of —1 or 1 means perfect 
correlation. A slight change in the same direction on the fitted curve can be achieved by either 
parameter. As an example we consider the high correlation of —0.95 between the hole mobility 
jip and the workfunction pc. This means that the current can be increased by either increasing 
the mobility up or lowering the workfunction pc. Such correlations do physically make sense, 
however they can complicate parameter fitting. Overall, our electrical parameter extraction 
example demonstrates, that rather complex physical models can be fitted to current-voltage 
curves thus making an efficient numerical analysis of measured data suitable. 


5. Conclusions 


We have given an overview on state-of-the-art electronic and optical device models for 
OLEDs. We demonstrated that numerical models are able to simulate different kinds 
of measurement setups for OLED characterization. We investigate the influence of the 
EGDM and the Gaussian density of states on the current-voltage curves, the dark-injection 
transients and the frequency response for a polymer device. All three characterization 
curves change significantly, if the disorder model ingredients are taken into account. The 
simulation results were successfully compared with analytical solutions, where applicable. 
We also show that the physical models can be combined with a nonlinear least-square 
fitting algorithm for extracting emission profiles and charge transport parameters from 
(hypothetical) measurements. In terms of optical device characteristics, quantitative 
outcoupling mode contributions and angular characteristics are presented as well as results 
of an emission zone extraction method. Simulations of electronic and optical processees in 
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a doped, multilayer OLED are also presented, thus highlighting, that comprehensive OLED 
device simulations are becoming feasible. 
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